Le monde ne suffit pas …

Une chaîne de traitement

Ce cours propose de décomposer les étapes successives d’une démarche de datamining en partant d’un exemple simple : la création d’une base de donnée sur les pays du Monde.
Claude Grasland (Professeur à l’Université de Paris)


Introduction

1 COLLECTE

1.1 L’API de la Banque Mondiale

Supposons que l’on souhaite télécharger la population, la superficie et le PIB des pays du monde. Plutôt que d’aller chercher des fichiers sur un site web, nous allons utiliser une API proposée par la Banque Mondial qui permet de télécharger les données facilement et surtout de les mettre à jour régulièrement. Pour cela on va installer le package R correspondant à l’API wbstats de la Banque mondiale.

https://cran.r-project.org/web/packages/wbstats/vignettes/Using_the_wbstats_package.html

Au moment du chargement du package, il est créé un fichier wb_cachelist qui fournit l’ensemble des donnes disponibles sous la forme d’une liste de tableaux de méta-données.

cat<-wb_cachelist
str(cat,max.level = 1)
List of 8
 $ countries    : tibble [304 × 18] (S3: tbl_df/tbl/data.frame)
 $ indicators   : tibble [16,649 × 8] (S3: tbl_df/tbl/data.frame)
 $ sources      : tibble [63 × 9] (S3: tbl_df/tbl/data.frame)
 $ topics       : tibble [21 × 3] (S3: tbl_df/tbl/data.frame)
 $ regions      : tibble [48 × 4] (S3: tbl_df/tbl/data.frame)
 $ income_levels: tibble [7 × 3] (S3: tbl_df/tbl/data.frame)
 $ lending_types: tibble [4 × 3] (S3: tbl_df/tbl/data.frame)
 $ languages    : tibble [23 × 3] (S3: tbl_df/tbl/data.frame)

1.1.1 Le tableau “countries”

Il fournit des renseignements de base sur les différents pays, leurs codes, etc.

str(cat$countries)
tibble [304 × 18] (S3: tbl_df/tbl/data.frame)
 $ iso3c             : chr [1:304] "ABW" "AFG" "AFR" "AGO" ...
 $ iso2c             : chr [1:304] "AW" "AF" "A9" "AO" ...
 $ country           : chr [1:304] "Aruba" "Afghanistan" "Africa" "Angola" ...
 $ capital_city      : chr [1:304] "Oranjestad" "Kabul" NA "Luanda" ...
 $ longitude         : num [1:304] -70 69.2 NA 13.2 19.8 ...
 $ latitude          : num [1:304] 12.52 34.52 NA -8.81 41.33 ...
 $ region_iso3c      : chr [1:304] "LCN" "SAS" NA "SSF" ...
 $ region_iso2c      : chr [1:304] "ZJ" "8S" NA "ZG" ...
 $ region            : chr [1:304] "Latin America & Caribbean" "South Asia" "Aggregates" "Sub-Saharan Africa" ...
 $ admin_region_iso3c: chr [1:304] NA "SAS" NA "SSA" ...
 $ admin_region_iso2c: chr [1:304] NA "8S" NA "ZF" ...
 $ admin_region      : chr [1:304] NA "South Asia" NA "Sub-Saharan Africa (excluding high income)" ...
 $ income_level_iso3c: chr [1:304] "HIC" "LIC" NA "LMC" ...
 $ income_level_iso2c: chr [1:304] "XD" "XM" NA "XN" ...
 $ income_level      : chr [1:304] "High income" "Low income" "Aggregates" "Lower middle income" ...
 $ lending_type_iso3c: chr [1:304] "LNX" "IDX" NA "IBD" ...
 $ lending_type_iso2c: chr [1:304] "XX" "XI" NA "XF" ...
 $ lending_type      : chr [1:304] "Not classified" "IDA" "Aggregates" "IBRD" ...

Le tableau comporte 304 observation et il mélange des pays (France), des fragments de pays (Réunion) et des agrégats de pays (Europe). Il faudra donc bien faire attention lors de l’extraction à réfléchir à ce que l’on souhaite utiliser. Par exemple, si l’on veut juste les pays :

pays<-cat$countries[cat$countries$income_level!="Aggregates",c("iso3c", "country","capital_city","longitude","latitude", "region","income_level")]
kable(head(pays))
iso3c country capital_city longitude latitude region income_level
ABW Aruba Oranjestad -70.0167 12.51670 Latin America & Caribbean High income
AFG Afghanistan Kabul 69.1761 34.52280 South Asia Low income
AGO Angola Luanda 13.2420 -8.81155 Sub-Saharan Africa Lower middle income
ALB Albania Tirane 19.8172 41.33170 Europe & Central Asia Upper middle income
AND Andorra Andorra la Vella 1.5218 42.50750 Europe & Central Asia High income
ARE United Arab Emirates Abu Dhabi 54.3705 24.47640 Middle East & North Africa High income

1.1.2 Le tableau indicators

Il comporta pas loin de 17000 variables … Autant dire qu’il est difficile de l’explorer facilement si l’on ne sait pas ce que l’on cherche.

indic<-cat$indicators
kable(head(indic,3))
indicator_id indicator unit indicator_desc source_org topics source_id source
1.0.HCount.1.90usd Poverty Headcount ($1.90 a day) NA The poverty headcount index measures the proportion of the population with daily per capita income (in 2011 PPP) below the poverty line. LAC Equity Lab tabulations of SEDLAC (CEDLAS and the World Bank). 11 , Poverty 37 LAC Equity Lab
1.0.HCount.2.5usd Poverty Headcount ($2.50 a day) NA The poverty headcount index measures the proportion of the population with daily per capita income (in 2005 PPP) below the poverty line. LAC Equity Lab tabulations of SEDLAC (CEDLAS and the World Bank). 11 , Poverty 37 LAC Equity Lab
1.0.HCount.Mid10to50 Middle Class ($10-50 a day) Headcount NA The poverty headcount index measures the proportion of the population with daily per capita income (in 2005 PPP) below the poverty line. LAC Equity Lab tabulations of SEDLAC (CEDLAS and the World Bank). 11 , Poverty 37 LAC Equity Lab

1.1.3 La recherche d’indicateurs

Supposons qu’on recherche les données récentes sur les émissions de CO2. On va utiliser le mot-clé CO2 pour rechercher les variables correspondantes dans le catalogue, ce qui donne 45 réponses

vars <- wbsearch(pattern = "CO2",fields="indicator")
kable(vars)
indicatorID indicator
5294 EN.ATM.CO2E.CP.KT CO2 emissions from cement production (thousand metric tons)
5295 EN.ATM.CO2E.EG.ZS CO2 intensity (kg per kg of oil equivalent energy use)
5296 EN.ATM.CO2E.FF.KT CO2 emissions from fossil-fuels, total (thousand metric tons)
5297 EN.ATM.CO2E.FF.ZS CO2 emissions from fossil-fuels (% of total)
5298 EN.ATM.CO2E.GDP CO2 emissions, industrial (kg per 1987 US$ of GDP)
5299 EN.ATM.CO2E.GF.KT CO2 emissions from gaseous fuel consumption (kt)
5300 EN.ATM.CO2E.GF.ZS CO2 emissions from gaseous fuel consumption (% of total)
5301 EN.ATM.CO2E.GL.KT CO2 emissions from gas flaring (thousand metric tons)
5302 EN.ATM.CO2E.KD.87.GD CO2 emissions, industrial (kg per 1987 US$ of GDP)
5303 EN.ATM.CO2E.KD.GD CO2 emissions (kg per 2010 US$ of GDP)
5304 EN.ATM.CO2E.KT CO2 emissions (kt)
5305 EN.ATM.CO2E.LF.KT CO2 emissions from liquid fuel consumption (kt)
5306 EN.ATM.CO2E.LF.ZS CO2 emissions from liquid fuel consumption (% of total)
5307 EN.ATM.CO2E.PC CO2 emissions (metric tons per capita)
5308 EN.ATM.CO2E.PP.GD CO2 emissions (kg per PPP $ of GDP)
5309 EN.ATM.CO2E.PP.GD.KD CO2 emissions (kg per 2017 PPP $ of GDP)
5310 EN.ATM.CO2E.SF.KT CO2 emissions from solid fuel consumption (kt)
5311 EN.ATM.CO2E.SF.ZS CO2 emissions from solid fuel consumption (% of total)
5312 EN.ATM.GHGO.KT.CE Other greenhouse gas emissions, HFC, PFC and SF6 (thousand metric tons of CO2 equivalent)
5314 EN.ATM.GHGT.KT.CE Total greenhouse gas emissions (kt of CO2 equivalent)
5316 EN.ATM.HFCG.KT.CE HFC gas emissions (thousand metric tons of CO2 equivalent)
5317 EN.ATM.METH.AG.KT.CE Agricultural methane emissions (thousand metric tons of CO2 equivalent)
5319 EN.ATM.METH.EG.KT.CE Methane emissions in energy sector (thousand metric tons of CO2 equivalent)
5322 EN.ATM.METH.KT.CE Methane emissions (kt of CO2 equivalent)
5323 EN.ATM.METH.PC Methane emissions (kt of CO2 equivalent per capita)
5325 EN.ATM.NOXE.AG.KT.CE Agricultural nitrous oxide emissions (thousand metric tons of CO2 equivalent)
5327 EN.ATM.NOXE.EG.KT.CE Nitrous oxide emissions in energy sector (thousand metric tons of CO2 equivalent)
5330 EN.ATM.NOXE.IN.KT.CE Industrial nitrous oxide emissions (thousand metric tons of CO2 equivalent)
5332 EN.ATM.NOXE.KT.CE Nitrous oxide emissions (thousand metric tons of CO2 equivalent)
5333 EN.ATM.NOXE.MT.CE Nitrous oxide emissions (metric tons of CO2 equivalent)
5334 EN.ATM.NOXE.PC Nitrous oxide emissions (metric tons of CO2 equivalent per capita)
5336 EN.ATM.PFCG.KT.CE PFC gas emissions (thousand metric tons of CO2 equivalent)
5343 EN.ATM.SF6G.KT.CE SF6 gas emissions (thousand metric tons of CO2 equivalent)
5348 EN.CLC.GHGR.MT.CE GHG net emissions/removals by LUCF (Mt of CO2 equivalent)
5353 EN.CO2.BLDG.MT CO2 emissions from residential buildings and commercial and public services (million metric tons)
5354 EN.CO2.BLDG.ZS CO2 emissions from residential buildings and commercial and public services (% of total fuel combustion)
5355 EN.CO2.ETOT.MT CO2 emissions from electricity and heat production, total (million metric tons)
5356 EN.CO2.ETOT.ZS CO2 emissions from electricity and heat production, total (% of total fuel combustion)
5357 EN.CO2.MANF.MT CO2 emissions from manufacturing industries and construction (million metric tons)
5358 EN.CO2.MANF.ZS CO2 emissions from manufacturing industries and construction (% of total fuel combustion)
5359 EN.CO2.OTHX.MT CO2 emissions from other sectors, excluding residential buildings and commercial and public services (million metric tons)
5360 EN.CO2.OTHX.ZS CO2 emissions from other sectors, excluding residential buildings and commercial and public services (% of total fuel combustion)
5361 EN.CO2.TRAN.MT CO2 emissions from transport (million metric tons)
5362 EN.CO2.TRAN.ZS CO2 emissions from transport (% of total fuel combustion)
8231 IN.ENV.CO2.CONC CO2 Emission (in thousand metric tons of Carbon)

On va finalement retenir deux indicateurs généraux

  • EN.ATM.CO2E.KT : émissions de CO2 en kilotonnes
  • EN.ATM.CO2E.PC : émissions de CO2 en tonnes par habitant

Puis examiner plus en détail leurs métadonnées

meta<-cat$indicators[cat$indicators$indicator_id %in% c("EN.ATM.CO2E.KT","EN.ATM.CO2E.PC"),]
kable(meta)
indicator_id indicator unit indicator_desc source_org topics source_id source
EN.ATM.CO2E.KT CO2 emissions (kt) NA Carbon dioxide emissions are those stemming from the burning of fossil fuels and the manufacture of cement. They include carbon dioxide produced during consumption of solid, liquid, and gas fuels and gas flaring. Carbon Dioxide Information Analysis Center, Environmental Sciences Division, Oak Ridge National Laboratory, Tennessee, United States. 19 , 6 , Climate Change, Environment 2 World Development Indicators
EN.ATM.CO2E.PC CO2 emissions (metric tons per capita) NA Carbon dioxide emissions are those stemming from the burning of fossil fuels and the manufacture of cement. They include carbon dioxide produced during consumption of solid, liquid, and gas fuels and gas flaring. Carbon Dioxide Information Analysis Center, Environmental Sciences Division, Oak Ridge National Laboratory, Tennessee, United States. 19 , 6 , Climate Change, Environment 2 World Development Indicators

1.1.4 L’extraction des données

Elle se fait à l’aide de la fonction wb_data qui comporte de nombreuses options. On peut par exemple extraire la valeur la plus récentes à l’aide de l’option mrv = 1

On peut ainsi établir le palmares des plus gros émetteurs de COE en valeur absolue …

tabCO2 <- wb_data(indicator = c("EN.ATM.CO2E.KT"),
                 return_wide = TRUE,
                 mrv=1,
                 country ="countries_only")



tabCO2<-tabCO2[order(tabCO2$EN.ATM.CO2E.KT,decreasing = T),]
kable(head(tabCO2,10))
iso2c iso3c country date EN.ATM.CO2E.KT unit obs_status footnote last_updated
CN CHN China 2016 9893038.0 NA NA NA 2020-12-16
US USA United States 2016 5006302.1 NA NA NA 2020-12-16
IN IND India 2016 2407671.5 NA NA NA 2020-12-16
RU RUS Russian Federation 2016 1732026.8 NA NA NA 2020-12-16
JP JPN Japan 2016 1135886.3 NA NA NA 2020-12-16
DE DEU Germany 2016 727972.8 NA NA NA 2020-12-16
IR IRN Iran, Islamic Rep. 2016 661710.2 NA NA NA 2020-12-16
KR KOR Korea, Rep. 2016 620302.4 NA NA NA 2020-12-16
SA SAU Saudi Arabia 2016 563449.2 NA NA NA 2020-12-16
ID IDN Indonesia 2016 563324.5 NA NA NA 2020-12-16

… ou bien par habitant.

tabCO2 <- wb_data(indicator = c("EN.ATM.CO2E.PC"),
                 return_wide = TRUE,
                 mrv=1,
                 country ="countries_only")



tabCO2<-tabCO2[order(tabCO2$EN.ATM.CO2E.PC,decreasing = T),]
kable(head(tabCO2,10))
iso2c iso3c country date EN.ATM.CO2E.PC unit obs_status footnote last_updated
QA QAT Qatar 2016 38.90147 NA NA NA 2020-12-16
CW CUW Curacao 2016 33.76146 NA NA NA 2020-12-16
TT TTO Trinidad and Tobago 2016 31.84485 NA NA NA 2020-12-16
KW KWT Kuwait 2016 24.95251 NA NA NA 2020-12-16
BH BHR Bahrain 2016 22.22898 NA NA NA 2020-12-16
AE ARE United Arab Emirates 2016 22.04083 NA NA NA 2020-12-16
NC NCL New Caledonia 2016 19.26650 NA NA NA 2020-12-16
GI GIB Gibraltar 2016 18.80401 NA NA NA 2020-12-16
BN BRN Brunei Darussalam 2016 18.25638 NA NA NA 2020-12-16
SA SAU Saudi Arabia 2016 17.36759 NA NA NA 2020-12-16

1.2 Puissance pays du Monde (1990-2020)

Nous allons essayer de constituer un tableau de la puissance des pays du Monde au cours de la période 1995-2020, à l’aide de six indicateurs de stock correspondant à différentes formes de puissance :

  • Puissance territoriale
    • SRF : Superficie totale du pays en km2
    • ARB : Superficie de terres arables en hectares
  • Puissance démographique
    • POP : Population totale en habitants
    • URB : Population urbaine en habitants
  • Puissance économique
    • GDP : Produit Intérieur Brut en parité de pouvoir d’achat
    • CO2 : Emissions de CO2 en tonnes

Ces indicateurs ont été choisis en raison de leur simplicité qui en assure la disponibilité pour la plupart des pays et pour la plupart des dates (excepté dans le cas du CO2 qui n’est mesuré que tardivement.)

1.2.1 Choix des indicateurs

On choisit un ensemble de données dont on connait l’identifiant et que l’on souhaite pouvoir analyser sur une période de temps longue.

world_data <- wb_data(indicator = c("AG.SRF.TOTL.K2","AG.LND.ARBL.HA","SP.POP.TOTL","SP.URB.TOTL","NY.GDP.MKTP.CD", "EN.ATM.CO2E.KT"),
                 return_wide = TRUE,
                 start_date = 1990,
                 end_date = 2019,
                 country ="countries_only")

world_data<-world_data[,-c(1,3)]

# recodage (attention : ordre alphabetique)
names(world_data)<-c("iso3c","date","ARB","SRF","CO2","GDP","POP","URB")

1.2.2 Recoller avec les données pays

On rajoute quelques données relatives au pays qui pourront être utiles par la suite.

pays<-cat$countries[cat$countries$income_level!="Aggregates",c("iso3c", "country","capital_city","longitude","latitude", "region","income_level")]
tab<-merge(pays,world_data, by="iso3c",all.x=F,all.y=T)
kable(head(tab))
iso3c country capital_city longitude latitude region income_level date ARB SRF CO2 GDP POP URB
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1990 2000 180 487.711 764887117 62149 31273
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1991 2000 180 531.715 872138715 64622 32507
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1992 2000 180 539.049 958463184 68235 34116
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1993 2000 180 649.059 1082979721 72504 35953
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1994 2000 180 660.060 1245688268 76700 37719
ABW Aruba Oranjestad -70.0167 12.5167 Latin America & Caribbean High income 1995 2000 180 707.731 1320474860 80324 39172
saveRDS(tab,"data/wb_don_1990_2019.Rdata")

1.2.3 Ajouter des métadonnées

# Extract meta
meta<-cat$indicators[cat$indicators$indicator_id %in% c("AG.SRF.TOTL.K2","AG.LND.ARBL.HA","SP.POP.TOTL","SP.URB.TOTL","NY.GDP.MKTP.CD", "EN.ATM.CO2E.KT"),]

# Select information and add personal code
meta<-data.frame(meta[,c(1,2,4,5)])
meta$shortcode<-c("ARB","SRF","CO2","GDP","POP","URB")
meta<-meta[,c(5,1,2,3,4)]

# Display
kable(meta)
shortcode indicator_id indicator indicator_desc source_org
ARB AG.LND.ARBL.HA Arable land (hectares) Arable land (in hectares) includes land defined by the FAO as land under temporary crops (double-cropped areas are counted once), temporary meadows for mowing or for pasture, land under market or kitchen gardens, and land temporarily fallow. Land abandoned as a result of shifting cultivation is excluded. Food and Agriculture Organization, electronic files and web site.
SRF AG.SRF.TOTL.K2 Surface area (sq. km) Surface area is a country’s total area, including areas under inland bodies of water and some coastal waterways. Food and Agriculture Organization, electronic files and web site.
CO2 EN.ATM.CO2E.KT CO2 emissions (kt) Carbon dioxide emissions are those stemming from the burning of fossil fuels and the manufacture of cement. They include carbon dioxide produced during consumption of solid, liquid, and gas fuels and gas flaring. Carbon Dioxide Information Analysis Center, Environmental Sciences Division, Oak Ridge National Laboratory, Tennessee, United States.
GDP NY.GDP.MKTP.CD GDP (current US$) GDP at purchaser’s prices is the sum of gross value added by all resident producers in the economy plus any product taxes and minus any subsidies not included in the value of the products. It is calculated without making deductions for depreciation of fabricated assets or for depletion and degradation of natural resources. Data are in current U.S. dollars. Dollar figures for GDP are converted from domestic currencies using single year official exchange rates. For a few countries where the official exchange rate does not reflect the rate effectively applied to actual foreign exchange transactions, an alternative conversion factor is used. World Bank national accounts data, and OECD National Accounts data files.
POP SP.POP.TOTL Population, total Total population is based on the de facto definition of population, which counts all residents regardless of legal status or citizenship. The values shown are midyear estimates. (1) United Nations Population Division. World Population Prospects: 2019 Revision. (2) Census reports and other statistical publications from national statistical offices, (3) Eurostat: Demographic Statistics, (4) United Nations Statistical Division. Population and Vital Statistics Reprot (various years), (5) U.S. Census Bureau: International Database, and (6) Secretariat of the Pacific Community: Statistics and Demography Programme.
URB SP.URB.TOTL Urban population Urban population refers to people living in urban areas as defined by national statistical offices. It is calculated using World Bank population estimates and urban ratios from the United Nations World Urbanization Prospects. Aggregation of urban and rural population may not add up to total population because of different country coverages. World Bank staff estimates based on the United Nations Population Division’s World Urbanization Prospects: 2018 Revision.
# Save meta
saveRDS(meta,"data/wb_don_1990_2019_meta.Rdata")
write.table(meta,"data/wb_don_1990_2019_meta.csv",sep=";",dec=",", row.names = F,fileEncoding = "UTF-8")

2 ESTIMATIONS

Le tableau que nous avons construit comporte encore de nombreuses valeurs manquantes. Or, notre objectif est de calculer la part d’un pays dans le total mondial ce qui n’est pas possibl si on ne dipose pas d’une estimation des valeurs de chacun des pays.

Nous allons donc construire un nouveau tableau où l’on essayera de remplir le maximum de valeurs manquantes tout en précisant la méthode d’estimatiopn utilisée.

2.1 Diagnostic des valeurs manquantes

Avant toute chose, nous allons estimer pour chacune des variables la part des valeurs manquantes en fonction des dates ou des pays.

2.1.1 Tableau de synthèse

On crée un tableau de synthèse des valeurs manquantes en mode colonne (toutes les variables les unes au dessus des autres et non pas côte à côte) et on utilise pour cela le package tidyverse qui est adapté à ce type de manipulation, notamment avec la fonction gather.

# Load data
don<-readRDS("data/wb_don_1990_2019.Rdata")



tabmis<-don %>% select(iso3c,country,date,ARB,SRF,CO2,GDP,POP,URB)  %>%
  gather("ARB", "SRF","CO2","GDP","POP","URB",key = VAR, value = DON) %>%
  mutate(Missing = is.na(DON))

kable(head(tabmis))             
iso3c country date VAR DON Missing
ABW Aruba 1990 ARB 2000 FALSE
ABW Aruba 1991 ARB 2000 FALSE
ABW Aruba 1992 ARB 2000 FALSE
ABW Aruba 1993 ARB 2000 FALSE
ABW Aruba 1994 ARB 2000 FALSE
ABW Aruba 1995 ARB 2000 FALSE

2.1.2 Valeurs manquantes par date

On réalise un graphique montrant le % de données manquantes par date et par variable en se servant de ggplot2

res <- tabmis %>% group_by(VAR, date) %>% 
                  summarise(nbmis=sum(Missing), nb=n(), pct = 100*nbmis/nb) 

p<-ggplot(res) + 
  aes(x=date,y=pct,color=VAR) +
  geom_line() +
  scale_y_continuous("% de pays à valeurs manquantes", breaks = c(0,10,20,30,40,50,60,70,80,90,100))+
  scale_x_continuous("Année", breaks = c(1990, 1995, 2000, 2005, 2010, 2015, 2020))
p

Les valeurs manquantes des années 1990 sont liées souvent à des pays qui n’existent pas encore comme le Sud-Soudan ou le Kosovo et qui ne disposent donc pas de données. Le fait qu’il demeure toujours 5% de pays non renseignés est lié à des pays de très petite taille mal documentés.

On observe que les données relatives au CO2 et aux terres arables (ARB) n’existent pour aucun pays après 2016. Il en va de même pour la superficie totale (SRF) des pays en 2019.

2.1.3 Valeurs manquantes par pays

On examine maintenant les pays ayant la plus forte proportion de valeurs manquantesn indépendamment de la variable concernée.

res <- tabmis %>% group_by(iso3c, country) %>% 
                  summarise(nbmis=sum(Missing), nb=n(), pct = 100*nbmis/nb) %>%
                  arrange(-pct)

kable(head(res,20))
iso3c country nbmis nb pct
MAF St. Martin (French part) 122 180 67.77778
SSD South Sudan 107 180 59.44444
SXM Sint Maarten (Dutch part) 97 180 53.88889
XKX Kosovo 92 180 51.11111
CUW Curacao 78 180 43.33333
GIB Gibraltar 64 180 35.55556
MCO Monaco 63 180 35.00000
SDN Sudan 57 180 31.66667
CHI Channel Islands 54 180 30.00000
NRU Nauru 54 180 30.00000
MNP Northern Mariana Islands 49 180 27.22222
MNE Montenegro 48 180 26.66667
VIR Virgin Islands (U.S.) 48 180 26.66667
ASM American Samoa 47 180 26.11111
GUM Guam 47 180 26.11111
SMR San Marino 44 180 24.44444
SRB Serbia 44 180 24.44444
ERI Eritrea 40 180 22.22222
IMN Isle of Man 40 180 22.22222
PRK Korea, Dem. People’s Rep. 37 180 20.55556

Trois cas apparaissent clairement :

  • Micro-états et territoires dépendants comme Monaco, Nauru, Guam, etc…

  • Etats issus de recompositions frontalières comme le Nord et le Sud-Soudan, le Kosovo, la Serbie, l’Erythrée … Ce sont par définition des pays qui n’existent pas à toutes les dates et l’attribution de valeurs dans le passé est une reconstitution.

  • Etats en crise ou dictatures qui ne fournissent pas de données comme la Corée du Nord ou sont dans l’incapacité de le faire comme la Syrie au cours des dernières années.

Si le premier cas n’est pas trop gênant (les petits états pèsent peu dans le total mondial), les second et troisième cas sont plus ennuyeux car il speuvent fausser le calcul de la part du total mondial des autres pays. On va donc tenter de proposer une estimation des valeurs manquantes qui permette de reconstituer le total mondial.

2.2 Estimation des valeurs manquantes

2.2.1 Fonction d’estimation

Nous avons construit ici une fonction d’estimation complexe qui utilise trois méthodes différentes selon la disposition des données manquantes :

  • interpolation : dans le cas d’une série interrompue sur un intervalle
  • extrapolation : dans le cas d’ue série où il manque les dernières valeurs et pour laquelle on applique la tendance moyenne des pays à valeurs non manquantes.
  • retropolation : dans le cas d’une série où il manque les premières valeurs et pour laquelle on applique la tendance moyenne des pays à valeurs non manquantes.

La méthode n’est toutefois pas applicable dans deux cas :

  1. absence complète de données pour un pays : car on ne dispose d’aucun point de calage
  2. absence complète de données pour une année : car on ne dispose pas de tendance moyenne de référence.
# load a data.frame with columns space, time, indic 
estim<-function(df   = don)               # dataframe
 
{


# Select variable
sel<-dcast(df,space~time)
M<-as.matrix(sel[,-1])
rownames(M)<-sel$space

####################################################
########   Estimation Procedure ####################
####################################################

# check for lines or column  filled with missing values and eliminate
check1<-apply(!is.na(M),FUN="sum",1)
check2<-apply(!is.na(M),FUN="sum",2)
M<-M[check1>0,check2>0]



dim(M)


# define dimensions
nr<-nrow(M)
nr
nc<-ncol(M)
nc


# create time matrix with NA
Mt<-matrix(rep(1:nc,nr),nrow=nr,nc=nc,byrow=T)
rownames(Mt)<-rownames(M)
colnames(Mt)<-colnames(M)
Mt[is.na(M)]<-NA
Mt[1:4,1:4]



# identify previous and next available valuer
Mt_min<-Mt
Mt_max<-Mt
for ( i in 1:nr) { 
  for ( j in 1:nc) { 
    if (is.na(Mt[i,j])) {
      Mt_min[i,j]<-max(Mt[i,1:j],na.rm=T)
      Mt_max[i,j]<-min(Mt[i,j:nc],na.rm=T)    
    }  
  }
}
Mt_min[1:4,1:4]
Mt_max[1:4,1:4]


# choose estimation method
M_met<-matrix("OK",nrow=nr,ncol=nc)
rownames(M_met)<-rownames(M)
colnames(M_met)<-colnames(M)
for ( i in 1:nr) { 
  for ( j in 1:nc) { 
    if (is.na(Mt[i,j])) {
      M_met[i,j]<-"IN"
      if (is.infinite(Mt_max[i,j])) {M_met[i,j]<-"EX"} 
      if (is.infinite(Mt_min[i,j])) {M_met[i,j]<-"RE"}  
    }   
  }
}
M_met[1:4,1:4]

###### Estimation of missing values  ####

M_est<-M
str(M_est)

#### step 1: Interpolation ######

for ( i in 1:nr) { 
  for ( j in 1:nc) { 
    if ((M_met[i,j]=="IN")) {
      t0<-Mt_min[i,j]
      t1<-Mt_max[i,j] 
      tacm<-(M[i,t1]/M[i,t0])**(1/(t1-t0))
      M_est[i,j]<-M[i,t0]*(tacm**(j-t0))
    }   
  }
}
M_est[1:4,1:4]


#### step 2: Extrapolation ######

for ( i in 1:nr) { 
  for ( j in 2:nc) { 
    if ((M_met[i,j]=="EX")) {
      N<-M_est
      t0<-j-1
      t1<-j
      N<-M_est[is.na(M_est[,t0])==F,]
      N<-N[is.na(N[,t1])==F,]
      tacm<-(sum(N[,t1])/sum(N[,t0]))**(1/(t1-t0)) 
      M_est[i,j]<-M_est[i,t0]*(tacm**(t1-t0))
    }   
  }
}
M_est[1:4,1:4]


#### step 3: Retropolation ######

for ( i in 1:nr) { 
  for ( j in (nc-1):1) { 
    if ((M_met[i,j]=="RE")) {
      N<-M_est
      t0<-j
      t1<-j+1
      N<-M_est[is.na(M_est[,t0])==F,]
      N<-N[is.na(N[,t1])==F,]
      tacm<-(sum(N[,t1])/sum(N[,t0]))**(1/(t1-t0)) 
      M_est[i,j]<-M_est[i,t1]/(tacm**(t1-t0))
    }   
  }
}
M_est[1:4,1:4]



#######################################################
############# EXPORT RESULTS  #########################
#######################################################

x<-reshape2::melt(M_est)
names(x)<-c("space","time","estim")
y<-reshape2::melt(M_met)
names(y)<-c("space","time","method")
z<-merge(x,y, by=c("space","time"))

return(z)
}

2.2.2 Application de la fonction

On va reconstituer pour chacune de nos variables les valeurs estimées lorsque cela est possible

# Load data
don<-readRDS("data/wb_don_1990_2019.Rdata")

# Estim ARB
df<-data.table(space=don$iso3c,time=don$date,index=don$ARB)
est<-estim(df)[,1:3]
 num [1:208, 1:27] 2000 7910000 2900000 579000 1000 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:208] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:27] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","ARB_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Estim SRF
df<-data.table(space=don$iso3c,time=don$date,index=don$SRF)
est<-estim(df)[,1:3]
 num [1:215, 1:29] 180 652860 1246700 28750 470 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:215] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:29] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","SRF_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Estim CO2
df<-data.table(space=don$iso3c,time=don$date,index=don$CO2)
est<-estim(df)[,1:3]
 num [1:207, 1:27] 488 2615 5115 5515 407 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:207] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:27] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","CO2_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Estim GDP
df<-data.table(space=don$iso3c,time=don$date,index=don$GDP)
est<-estim(df)[,1:3]
 num [1:213, 1:30] 764887117 NA 11228764963 2028553750 1029048482 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:213] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:30] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","GDP_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Estim POP
df<-data.table(space=don$iso3c,time=don$date,index=don$POP)
est<-estim(df)[,1:3]
 num [1:217, 1:30] 62149 12412308 11848386 3286542 54509 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:217] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:30] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","POP_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Estim URB
df<-data.table(space=don$iso3c,time=don$date,index=don$URB)
est<-estim(df)[,1:3]
 num [1:215, 1:30] 31273 2628554 4400964 1197222 51627 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:215] "ABW" "AFG" "AGO" "ALB" ...
  ..$ : chr [1:30] "1990" "1991" "1992" "1993" ...
names(est)<-c("iso3c","date","URB_est")
don<-merge(don,est, by=c("iso3c","date"), all.x=T,al.y=F)

# Save
saveRDS(don,"data/wb_don_1990_2019.Rdata")

Attention ! les valeurs estimées sont parfois très éloignées de la réalité, surtout dans le cas des extrapolations où elles suivent la tendance mondiale. Mais cette méthode d’estimation permet, comme nous l’avons expliqué, de pouvir calculer le total mondial et, du coup, de pouvoir estimer pour chaque critère la part du pays dans la population mondiale et son rang à chacune des dates.

3 TRANSFORMATION

On se propose d’enrichir notre tableau à l’aide de nouvelles variables

  1. Rang des pays pour chacun des critères
  2. Part des pays dans le total mondial pour chacun des critères
  3. Indicateurs synthétiques de puissance combinant les différents critères

3.1 Rangs mondial des pays

3.1.1 Calcul des rangs par variable et date

On utilise dplyr pour calculer le rand de chaque variable à chacune des dates. On ajoute une variable supplémenaire qui est la moyenne des rang sur les six classements.

don<-readRDS("data/wb_don_1990_2019.Rdata")
#don<-as_tibble(don)
don2<-don %>%  group_by(date) %>% mutate(ARB_rnk = rank(desc(ARB_est)),
                                         SRF_rnk = rank(desc(SRF_est)),
                                         CO2_rnk = rank(desc(CO2_est)),
                                         GDP_rnk = rank(desc(GDP_est)),
                                         POP_rnk = rank(desc(POP_est)),
                                         URB_rnk = rank(desc(URB_est)),
                                         TOT_rnk = (ARB_rnk+SRF_rnk+CO2_rnk+GDP_rnk+POP_rnk+URB_rnk)/6) 

3.1.2 Application au calcul des 10 pays les plus puissants en 1996 et 2016

tab<-don2 %>% filter(date==1996) %>% 
             select(iso3c, country, TOT_rnk, ARB_rnk,SRF_rnk,CO2_rnk,GDP_rnk,POP_rnk,URB_rnk,TOT_rnk) %>%
             arrange(TOT_rnk)
kable(head(tab,10), digits=1)
date iso3c country TOT_rnk ARB_rnk SRF_rnk CO2_rnk GDP_rnk POP_rnk URB_rnk
1996 USA United States 2.0 1 3 1 1 3 3
1996 CHN China 3.2 4 4 2 7 1 1
1996 IND India 5.7 2 7 6 15 2 2
1996 RUS Russian Federation 5.7 3 1 3 16 6 5
1996 BRA Brazil 7.3 5 5 17 8 5 4
1996 MEX Mexico 12.3 14 13 15 13 11 8
1996 IDN Indonesia 14.0 17 14 20 22 4 7
1996 CAN Canada 14.3 6 2 8 10 33 27
1996 FRA France 17.8 16 46 12 4 18 11
1996 DEU Germany 19.0 24 61 5 3 12 9
tab<-don2 %>% filter(date==2016) %>% 
             select(iso3c, country, TOT_rnk, ARB_rnk,SRF_rnk,CO2_rnk,GDP_rnk,POP_rnk,URB_rnk,TOT_rnk) %>%
             arrange(TOT_rnk)
kable(head(tab,10), digits=1)
date iso3c country TOT_rnk ARB_rnk SRF_rnk CO2_rnk GDP_rnk POP_rnk URB_rnk
2016 CHN China 2.2 4 4 1 2 1 1
2016 USA United States 2.3 2 3 2 1 3 3
2016 IND India 3.7 1 7 3 7 2 2
2016 RUS Russian Federation 6.0 3 1 4 12 9 7
2016 BRA Brazil 7.0 5 5 14 9 5 4
2016 IDN Indonesia 10.3 13 14 10 16 4 5
2016 MEX Mexico 12.2 14 13 12 15 11 8
2016 CAN Canada 16.2 7 2 11 10 38 29
2016 IRN Iran, Islamic Rep. 16.8 21 16 7 26 18 13
2016 TUR Turkey 18.8 15 35 17 17 17 12

3.2 Part du total mondial

Le poids d’un pays dans le monde est entendu ici comme sa part (son pourcentage) du total mondial. il correspond donc à une mesure de puissance relative.

3.2.1 Calcul des poids par variable et date

On utilise à nouveau dplyr pour calculer la part de chaque pays pour chaque variable à chacune des dates. On ajoute une variable supplémenaire qui est la moyenne des poids sur les six classements.

don2<-don2 %>%  group_by(date) %>% mutate(ARB_wgt = 100*ARB_est/sum(ARB_est,na.rm=T),
                                         SRF_wgt = 100*SRF_est/sum(SRF_est,na.rm=T),
                                         CO2_wgt = 100*CO2_est/sum(CO2_est,na.rm=T),
                                         GDP_wgt = 100*GDP_est/sum(GDP_est,na.rm=T),
                                         POP_wgt = 100*POP_est/sum(POP_est,na.rm=T),
                                         URB_wgt = 100*URB_est/sum(URB_est,na.rm=T),
                                         TOT_wgt = (ARB_wgt+SRF_wgt+CO2_wgt+GDP_wgt+POP_wgt+URB_wgt)/6)

3.2.2 Application au calcul des 10 pays les plus puissants en 1996 et 2016

tab<-don2 %>% filter(date==1996) %>% 
             select(iso3c, country, TOT_wgt, ARB_wgt,SRF_wgt,CO2_wgt,GDP_wgt,POP_wgt,URB_wgt,TOT_wgt) %>%
             arrange(desc(TOT_wgt))
kable(head(tab,10), digits=1)
date iso3c country TOT_wgt ARB_wgt SRF_wgt CO2_wgt GDP_wgt POP_wgt URB_wgt
1996 USA United States 13.6 12.8 7.3 23.1 25.8 4.7 8.0
1996 CHN China 11.6 8.5 7.3 15.3 2.8 21.1 14.9
1996 IND India 7.7 11.5 2.5 3.9 1.3 17.0 10.1
1996 RUS Russian Federation 6.2 9.0 13.0 7.1 1.3 2.6 4.2
1996 JPN Japan 4.6 0.3 0.3 5.3 15.4 2.2 3.8
1996 BRA Brazil 3.7 4.1 6.5 1.3 2.7 2.9 5.0
1996 DEU Germany 2.8 0.8 0.3 3.9 8.0 1.4 2.3
1996 CAN Canada 2.7 3.3 7.6 2.1 2.0 0.5 0.9
1996 AUS Australia 2.0 2.6 5.9 1.3 1.3 0.3 0.6
1996 FRA France 1.9 1.3 0.4 1.7 5.1 1.0 1.7
tab<-don2 %>% filter(date==2016) %>% 
             select(iso3c, country, TOT_wgt, ARB_wgt,SRF_wgt,CO2_wgt,GDP_wgt,POP_wgt,URB_wgt,TOT_wgt) %>%
             arrange(desc(TOT_wgt))
kable(head(tab,10), digits=1)
date iso3c country TOT_wgt ARB_wgt SRF_wgt CO2_wgt GDP_wgt POP_wgt URB_wgt
2016 CHN China 16.3 8.4 7.2 29.5 14.8 18.6 19.5
2016 USA United States 11.5 10.7 7.4 14.9 24.7 4.4 6.6
2016 IND India 8.8 11.0 2.5 7.2 3.0 17.9 10.9
2016 RUS Russian Federation 5.5 8.6 13.0 5.2 1.7 2.0 2.7
2016 BRA Brazil 3.8 5.7 6.5 1.4 2.4 2.8 4.4
2016 CAN Canada 2.6 3.1 7.6 1.6 2.0 0.5 0.7
2016 JPN Japan 2.5 0.3 0.3 3.4 6.5 1.7 2.9
2016 IDN Indonesia 2.2 1.7 1.4 1.7 1.2 3.5 3.5
2016 AUS Australia 2.1 3.2 5.9 1.1 1.6 0.3 0.5
2016 DEU Germany 1.8 0.8 0.3 2.2 4.6 1.1 1.6
saveRDS(don2,"data/wb_don_1990_2019.Rdata")
write.table(don2,"data/wb_don_1990_2019.csv",sep=";",dec=",", row.names = F,fileEncoding = "UTF-8")

4 SPATIALISATION

Au cours de cette étape, on va essayer d’associer à chaque pays du Monde une géométrie permettant d’en faire la cartographie. Longtemps compliquée, cette opération est maintenant facilitée par le packge sf (spatial features) qui permet grosso modo de stocker la géométrie (contour des pays) sous la forme d’une simple colonne ajoutée au tableau de données. On peut ensuite facilement réaliser des changements de projections et des cartes statiques ou dynamiques.

Une difficulté plus importante est de trouver un fonds de carte où le code des unités géométriques corresponde à celui des unités statistiques que nous avons collectées. Cette opération de jointure s’avère toujours délicate et elle l’est tout particulièrement dans le cas des pays du monde qui sont un objet finalement mal défini, tant sur le plan politique que sur le plan statistique.

4.1 L’API Natural Earth

Nous allons ici utiliser le fonds de carte Natural Earth qui est un fonds de carte libre de droit et mis à jour régulièrement. Le site web du projet se situe à l’adresse suivante :

https://www.naturalearthdata.com/

Il indique ses objectifs comme suit :

“Natural Earth is a public domain map dataset available at 1:10m, 1:50m, and 1:110 million scales. Featuring tightly integrated vector and raster data, with Natural Earth you can make a variety of visually pleasing, well-crafted maps with cartography or GIS software.[…] Natural Earth was built through a collaboration of many volunteers and is supported by NACIS (North American Cartographic Information Society), and is free for use in any type of project (see our Terms of Use page for more information).”

On peut télécharger les différents fonds de carte sur le site web, mais dans une perspective de mise à jour automatique régulière du fonds de carte il est plus pertinent d’utiliser l’API rnaturalearthqui permet d’accéder directement à la plupart des fonds de carte avec juste quelques lignes de code.

4.1.1 le fonds de carte countries110 (175 unités)

On va télécharger tout d’abord le fonds de carte des pays du Monde avec une forte généralisation des contours countries110 et le transformer en objet de type sf avant de le visualiser et d’ examiner le nombre d’unités

map<-st_as_sf(sovereignty110)
par(mar=c(0,0,0,0), mfrow=c(1,1))
plot(map$geometry,col="lightyellow")

dim(map)
[1] 171  64

Ce fonds de carte comporte 175 unités spatiales, mais de quoi s’agit-il exactement. Les métadonnées associées permettent de se faire une idée plus précise de la nature exacte de ces unités. Prenons pour cela quelques exempes

sel<-map[map$adm0_a3 %in% c("FRA", "NCL","ATA","ATF","USA", "PRI","CHN","TWN","MAR", "SAH","CHN","TWN","ISR","PSX"),c("sovereignt","sov_a3","type","admin", "adm0_a3","name","note_adm0","iso_a3","wb_a3")]
kable(sel)
sovereignt sov_a3 type admin adm0_a3 name note_adm0 iso_a3 wb_a3 geometry
6 Antarctica ATA Indeterminate Antarctica ATA Antarctica NA ATA NA MULTIPOLYGON (((-59.57209 -…
75 Israel ISR Sovereign country Israel ISR Israel NA ISR ISR MULTIPOLYGON (((35.8211 33….
96 Morocco MAR Sovereign country Morocco MAR Morocco NA MAR MAR MULTIPOLYGON (((-5.193863 3…
131 Western Sahara SAH Indeterminate Western Sahara SAH W. Sahara Self admin. ESH NA MULTIPOLYGON (((-8.794884 2…
157 Taiwan TWN Sovereign country Taiwan TWN Taiwan NA TWN NA MULTIPOLYGON (((121.7778 24…

Les exemples présentés dans le tableau ci-dessus montrent la complexité du problème de définition et de représentation cartographique des “pays” ou “bouts du monde”. Quelques remarques :

  1. La France (FR1) en tant qu’état souverain regroupe ici cartographiquement la partie métropolitaine du pays et les Départements d’Outre-Mer (Guyane Française, Réunion, Martinique, Guadeloupe) en une seule entité spatiale, mais elle met à part la Nouvelle Calédonie et les îles antarctiques.
  2. Porto Rico (PRI) est considéré comme une dépendance des Etats-Unis (US1) au même titre que la Nouvelle Calédonie(NCL) est considérée comme une dépendance de la France (FR1).
  3. Le Sahara occidental (SAH) est considéré comme une zone indéterminée bien qu’il soit occupé par le Maroc (MAR).
  4. la Palestine (PSX) est considéré comme une zone disputée mais rattachée en terme de souveraineté à Israël (ISR) et une note précise qu’elle est partiellement semi-administrée. Le code sur trois caractères des territoires palestiniens est très variable selon les organisations (PSX, PSE, WBG).
  5. Taïwan (TWN) est présenté comme un état souverain, mais son code ISO3 est manquant pour la banque mondiale car la Chine refuse de le reconnaître.
  6. Plusieurs états souverains de petite taille sont absents de ce fonds de carte qui ne regroupe que 175 unités soit moins que les 193 pays membres des Nations-Unies. La plupart des îles du Pacifique sont en particulier éliminées car leur surface les rendrait invisible pour le degré de généralisation cartographique adopté.

4.1.2 le fonds de carte sovereignty110 (171 unités)

On peut obtenir un fonds différent en installant le package complémentaire rnaturalearthdata qui permet notamment de distinguer le fonds de carte des countries (c’est-à-dire des “bouts du monde” souverains ou non) et des sovereignty (c’est-à-dire des états souverains)

map<-st_as_sf(sovereignty110)
par(mar=c(0,0,0,0))
plot(map$geometry,col="lightyellow")

dim(map)
[1] 171  64

Le fonds de carte permet désormais de récupérer la plupart des pays souverains du Monde, y compris les petits états insulaires du Pacifique, mais il fait disparaître de façon sélective les territoires indéterminés ou disputés. Ainsi, le Sahra Occidental demeure partiellement séparé du Maroc mais les territoires palestiniens sont annexés à Israël ainsi que le plateau du Golan ce qui n’est évidemment pas un choix neutred’un point de vue géoolitique.

par(mfrow=c(1,2))
plot(map[map$sov_a3 %in% c("ISR","JOR","SYR","LBN","EGY"),]$geometry, col=c("gray80","orange","gray80","gray80","gray80"))
title("Limits of Israël",cex=0.5)
plot(map[map$sov_a3 %in% c("MAR","SAH","DZA","MRT"),]$geometry, col=c("gray70","orange","gray70","lightyellow"))
title("Limits of Morocco")

4.1.3 Le fonds de carte tinycountries110

On peut aussi revenir au fonds de carte des countries et extraire les “petits pays” en ne conservant que leur point central, sans tracer un polygône de contour. On pourra ainsi les cartographier sous forme ponctuelle.

map<-st_as_sf(countries110)
small<-st_as_sf(tiny_countries110)
par(mar=c(0,0,0,0), mfrow=c(1,1))
plot(map$geometry,col="lightyellow")
plot(small$geometry,col="red", add=T)

4.1.4 Le fonds de carte countries50

On peut également choisir un fonds moins généralisé dans leque tous les petits pays seront présents

map<-st_as_sf(countries50)

par(mar=c(0,0,0,0))
plot(map$geometry,col="lightyellow")

4.1.5 Autres fonds de carte :

Il existe toute une série d’autres fonds de carte dans le package Natural Earth, notamment avec des résolutions plus précises, mais on se limitera ici à l’exploration des fonds de cart utile pour produire des cartes à contour généralisé couvrant le monde entier.

4.2 Application

Nous allons essayer de construire un fonds de carte qui permette de visualiser l’ensemble des données présentes dans le fichier de la banque mondiale en 2015. Plus précisément, nous allons construire deux fonds de carte, l’un avec une résolution faible ne comportant que 175 pays et l’autre avec une résolution détaillée comportant tous les pays.

4.2.1 Fonds de carte wb_map_low

On se limite aux plus grands pays

# Load map
map<-st_as_sf(countries110)
map<-map[c("adm0_a3","name")]
names(map)<-c("iso3c","name","geometry")

# Add polygons center
coo<-st_coordinates(st_centroid(map,of_largest_polygon = T))
map$Lon<-coo[,1]
map$Lat<-coo[,2]


# adjust some codes
map$iso3c[map$iso3c=="KOS"]<-"XKX"    # Kosovo
map$iso3c[map$iso3c=="PSX"]<-"PSE"    # Palestinian territories
map$iso3c[map$iso3c=="SDS"]<-"SSD"    # South Sudan

# Save
st_write(map,"data/wb_map_low.shp",delete_dsn=T)
Deleting source `data/wb_map_low.shp' failed
Writing layer `wb_map_low' to data source `data/wb_map_low.shp' using driver `ESRI Shapefile'
Writing 177 features with 4 fields and geometry type Multi Polygon.
saveRDS(map,"data/wb_map_low.Rdata")

4.2.2 Fonds de carte wb_map_high

# Load map
map<-st_as_sf(countries50)
map<-map[c("adm0_a3","name")]
names(map)<-c("iso3c","name","geometry")

# Add polygons center
coo<-st_coordinates(st_centroid(map,of_largest_polygon = T))
map$Lon<-coo[,1]
map$Lat<-coo[,2]


# adjust some codes
map$iso3c[map$iso3c=="KOS"]<-"XKX"    # Kosovo
map$iso3c[map$iso3c=="PSX"]<-"PSE"    # Palestinian territories
map$iso3c[map$iso3c=="SDS"]<-"SSD"    # South Sudan

# Save
st_write(map,"data/wb_map_high.shp",delete_dsn=T)
Deleting source `data/wb_map_high.shp' failed
Writing layer `wb_map_high' to data source `data/wb_map_high.shp' using driver `ESRI Shapefile'
Writing 241 features with 4 fields and geometry type Multi Polygon.
saveRDS(map,"data/wb_map_high.Rdata")

5 VISUALISATIONS

On va élaborer une série de représentations interactives qui pourront ultérieurement être reprises dans une application de type dashboard. On va utiliser pour cela principalement le package plotly.

L’exemple retenu sera celui des émissions de CO2 de 1995 à 2015

5.1 Courbes interactives

5.1.1 CO2

don<-readRDS("data/wb_don_1990_2019.Rdata")
tab<-don %>% filter(date > 1994, date < 2016) %>%
  select (iso3c, country, date, CO2_est, POP_est, CO2_rnk, CO2_wgt) %>% 
  filter(CO2_rnk <11) 

tab<-as.data.frame(tab)


p <- plot_ly(tab,
             x = ~date,
             y = ~CO2_wgt,
             color= ~country) %>%
             add_lines()%>%
         layout(title = "Top 10 des émissions de CO2 (1995-2015)",
         yaxis = list(title = "% des émissions mondiales de CO2", type="log"), 
         xaxis = list(title = "Année"))
p

5.1.2 PIB

don<-readRDS("data/wb_don_1990_2019.Rdata")
tab<-don %>% filter(date > 1994, date < 2016) %>%
  select (iso3c, country, date, GDP_est, POP_est, GDP_rnk, GDP_wgt) %>% 
  filter(GDP_rnk <21) 

tab<-as.data.frame(tab)


p <- plot_ly(tab,
             x = ~date,
             y = ~GDP_wgt,
             color= ~country) %>%
             add_lines()%>%
         layout(title = "Top 20 du PIB en ppa (1995-2015)",
         yaxis = list(title = "% du PIB mondial en ppa", type="log"), 
         xaxis = list(title = "Année"))
p

5.2 Cartes interactives

5.2.1 CO2

# Load map
map<-readRDS("data/wb_map_high.Rdata")

# Load don for 2015
don<-readRDS("data/wb_don_1990_2019.Rdata")
don<-don[don$date==2015,]

# Merge with map
mapdon<-merge(don,map,by="iso3c",all.x=F,all.y=T)

# Create map



g <- list(showframe = TRUE,
          framecolor= toRGB("gray50"),
          coastlinecolor = toRGB("black"),
          showland = TRUE,
          landcolor = toRGB("lightyellow"),
          showcountries = TRUE,
          countrycolor = toRGB("black"),
          countrywidth = 0.2,
 #        projection = list(type = 'azimuthal equal area'))
          projection = list(type = 'Mercator'))



p<- plot_geo(mapdon)%>%
  add_markers(x = ~Lon,
              y = ~Lat,
              sizes = c(0, 1000),
              size = ~CO2_wgt,
              color= "red",
              hoverinfo = "text",
              text = ~paste('Pays: ', country,
                            '<br /> % mondial : ',round(CO2_wgt,3),
                            '<br /> rang mondial : ',round(CO2_rnk,0))) %>%
  layout(geo = g,
         title = "Part des émissions mondiales de CO2 en 2015",
         width=800,
         height = 400)
p

5.2.2 GDP

# Load map
map<-readRDS("data/wb_map_high.Rdata")

# Load don for 2015
don<-readRDS("data/wb_don_1990_2019.Rdata")
don<-don[don$date==2015,]

# Merge with map
mapdon<-merge(don,map,by="iso3c",all.x=F,all.y=T)

# Create map



g <- list(showframe = TRUE,
          framecolor= toRGB("gray50"),
          coastlinecolor = toRGB("black"),
          showland = TRUE,
          landcolor = toRGB("lightyellow"),
          showcountries = TRUE,
          countrycolor = toRGB("black"),
          countrywidth = 0.2,
 #        projection = list(type = 'azimuthal equal area'))
          projection = list(type = 'Mercator'))



p<- plot_geo(mapdon)%>%
  add_markers(x = ~Lon,
              y = ~Lat,
              sizes = c(0, 1000),
              size = ~GDP_wgt,
              color= "red",
              hoverinfo = "text",
              text = ~paste('Pays: ', country,
                            '<br /> % mondial : ',round(GDP_wgt,3),
                            '<br /> rang mondial : ',round(GDP_rnk,0))) %>%
  layout(geo = g,
         title = "Part du PIB mondial (ppa) en 2015",
         width=800,
         height = 400)
p